{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [],
   "source": [
    "%matplotlib inline\n",
    "import numpy as np\n",
    "import pylab as pl\n",
    "from scipy import sparse\n",
    "from scipy.sparse import csgraph"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 稀疏矩阵-sparse"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 稀疏矩阵的储存形式"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "dict_keys([(2, 3), (3, 3), (4, 3)])\n",
      "dict_values([1.0, 2.0, 3.0])\n"
     ]
    }
   ],
   "source": [
    "from scipy import sparse\n",
    "a = sparse.dok_matrix((10, 5))\n",
    "\n",
    "a[2, 3] = 1.0\n",
    "a[3, 3] = 2.0\n",
    "a[4, 3] = 3.0\n",
    "print(a.keys())\n",
    "print(a.values())"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[list([]) list([]) list([1.0]) list([3.0, 2.0]) list([]) list([]) list([])\n",
      " list([]) list([]) list([])]\n",
      "[list([]) list([]) list([3]) list([2, 4]) list([]) list([]) list([])\n",
      " list([]) list([]) list([])]\n"
     ]
    }
   ],
   "source": [
    "b = sparse.lil_matrix((10, 5))\n",
    "b[2, 3] = 1.0\n",
    "b[3, 4] = 2.0\n",
    "b[3, 2] = 3.0\n",
    "print(b.data)\n",
    "print(b.rows)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[3 4 2 3] [2 3 3 2] [ 1  2  3 10]\n",
      "[[ 0  0  0  0  0  0]\n",
      " [ 0  0  0  0  0  0]\n",
      " [ 0  0  0 11  0  0]\n",
      " [ 0  0  3  0  2  0]\n",
      " [ 0  0  0  0  0  0]]\n"
     ]
    }
   ],
   "source": [
    "row = [2, 3, 3, 2]\n",
    "col = [3, 4, 2, 3]\n",
    "data = [1, 2, 3, 10]\n",
    "c = sparse.coo_matrix((data, (row, col)), shape=(5, 6))\n",
    "print (c.col, c.row, c.data)\n",
    "print (c.toarray())"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 最短路径"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {},
   "outputs": [
    {
     "ename": "SyntaxError",
     "evalue": "invalid syntax (<ipython-input-14-adc4f84eb492>, line 3)",
     "output_type": "error",
     "traceback": [
      "\u001b[1;36m  File \u001b[1;32m\"<ipython-input-14-adc4f84eb492>\"\u001b[1;36m, line \u001b[1;32m3\u001b[0m\n\u001b[1;33m    digraph graph1 {\u001b[0m\n\u001b[1;37m                 ^\u001b[0m\n\u001b[1;31mSyntaxError\u001b[0m\u001b[1;31m:\u001b[0m invalid syntax\n"
     ]
    }
   ],
   "source": [
    "\n",
    "#%figonly=使用稀疏矩阵可以表示有向图\n",
    "digraph graph1 {\n",
    "    rankdir=LR;\n",
    "    size=\"8,5\"\n",
    "    node [shape = circle];\n",
    "    A -> B [ label = \"10\" ];\n",
    "    B -> C [ label = \"5\" ];\n",
    "    A -> C [ label = \"3\" ];\n",
    "    C -> D [ label = \"7\" ];\n",
    "    D -> A [ label = \"4\" ];\n",
    "    D -> C [ label = \"6\" ];\n",
    "}"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "matrix([[ 0., 10.,  3.,  0.],\n",
       "        [ 0.,  0.,  5.,  0.],\n",
       "        [ 0.,  0.,  0.,  7.],\n",
       "        [ 4.,  0.,  6.,  0.]])"
      ]
     },
     "execution_count": 15,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "w = sparse.dok_matrix((4, 4))\n",
    "\n",
    "edges = [(0, 1, 10), (1, 2, 5), (0, 2, 3),\n",
    "         (2, 3,  7), (3, 0, 4), (3, 2, 6)]\n",
    "\n",
    "for i, j, v in edges:\n",
    "    w[i, j] = v\n",
    "\n",
    "w.todense()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "metadata": {},
   "outputs": [
    {
     "ename": "IndexError",
     "evalue": "index 0 is out of bounds for axis 0 with size 0",
     "output_type": "error",
     "traceback": [
      "\u001b[1;31m---------------------------------------------------------------------------\u001b[0m",
      "\u001b[1;31mIndexError\u001b[0m                                Traceback (most recent call last)",
      "\u001b[1;32m<ipython-input-17-c7c2ac13b02d>\u001b[0m in \u001b[0;36m<module>\u001b[1;34m()\u001b[0m\n\u001b[0;32m      5\u001b[0m \u001b[0mH\u001b[0m\u001b[1;33m,\u001b[0m \u001b[0mW\u001b[0m \u001b[1;33m=\u001b[0m \u001b[0mbimg\u001b[0m\u001b[1;33m.\u001b[0m\u001b[0mshape\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0;32m      6\u001b[0m \u001b[1;33m\u001b[0m\u001b[0m\n\u001b[1;32m----> 7\u001b[1;33m \u001b[0mx0\u001b[0m\u001b[1;33m,\u001b[0m \u001b[0mx1\u001b[0m \u001b[1;33m=\u001b[0m \u001b[0mnp\u001b[0m\u001b[1;33m.\u001b[0m\u001b[0mwhere\u001b[0m\u001b[1;33m(\u001b[0m\u001b[0mbimg\u001b[0m\u001b[1;33m[\u001b[0m\u001b[0mH\u001b[0m\u001b[1;33m//\u001b[0m\u001b[1;36m2\u001b[0m\u001b[1;33m,\u001b[0m \u001b[1;33m:\u001b[0m\u001b[1;33m]\u001b[0m\u001b[1;33m==\u001b[0m\u001b[1;36m0\u001b[0m\u001b[1;33m)\u001b[0m\u001b[1;33m[\u001b[0m\u001b[1;36m0\u001b[0m\u001b[1;33m]\u001b[0m\u001b[1;33m[\u001b[0m\u001b[1;33m[\u001b[0m\u001b[1;36m0\u001b[0m\u001b[1;33m,\u001b[0m \u001b[1;33m-\u001b[0m\u001b[1;36m1\u001b[0m\u001b[1;33m]\u001b[0m\u001b[1;33m]\u001b[0m \u001b[1;31m#❷\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0m\u001b[0;32m      8\u001b[0m \u001b[0mbimg\u001b[0m\u001b[1;33m[\u001b[0m\u001b[0mH\u001b[0m\u001b[1;33m//\u001b[0m\u001b[1;36m2\u001b[0m\u001b[1;33m,\u001b[0m \u001b[1;33m:\u001b[0m\u001b[0mx0\u001b[0m\u001b[1;33m]\u001b[0m \u001b[1;33m=\u001b[0m \u001b[1;36m0\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0;32m      9\u001b[0m \u001b[0mbimg\u001b[0m\u001b[1;33m[\u001b[0m\u001b[0mH\u001b[0m\u001b[1;33m//\u001b[0m\u001b[1;36m2\u001b[0m\u001b[1;33m,\u001b[0m \u001b[0mx1\u001b[0m\u001b[1;33m:\u001b[0m\u001b[1;33m]\u001b[0m \u001b[1;33m=\u001b[0m \u001b[1;36m0\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n",
      "\u001b[1;31mIndexError\u001b[0m: index 0 is out of bounds for axis 0 with size 0"
     ]
    }
   ],
   "source": [
    "img = pl.imread(\"maze.jpg\")\n",
    "sx, sy = (400, 979)\n",
    "ex, ey = (398,  25)\n",
    "bimg = np.all(img > 0.81, axis=2) #❶\n",
    "H, W = bimg.shape\n",
    "\n",
    "x0, x1 = np.where(bimg[H//2, :]==0)[0][[0, -1]] #❷\n",
    "bimg[H//2, :x0] = 0\n",
    "bimg[H//2, x1:] = 0"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "metadata": {},
   "outputs": [],
   "source": [
    "#上下相邻白色像素\n",
    "mask = (bimg[1:, :] & bimg[:-1, :]) \n",
    "idx = np.where(mask.ravel())[0]\n",
    "vedge = np.c_[idx, idx + W]\n",
    "pl.imsave(\"tmp.png\", mask, cmap=\"gray\")\n",
    "\n",
    "#左右相邻白色像素\n",
    "mask = (bimg[:, 1:] & bimg[:, :-1])\n",
    "y, x = np.where(mask)\n",
    "idx = y * W + x\n",
    "hedge = np.c_[idx, idx + 1]\n",
    "\n",
    "edges = np.vstack([vedge, hedge]) #❶\n",
    "\n",
    "values = np.ones(edges.shape[0])\n",
    "w = sparse.coo_matrix((values, (edges[:, 0], edges[:, 1])),  #❷\n",
    "                      shape=(bimg.size, bimg.size))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 19,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "(1, 801600)\n",
      "(1, 801600)\n"
     ]
    }
   ],
   "source": [
    "from scipy.sparse import csgraph\n",
    "startid = sy * W + sx\n",
    "endid   = ey * W + ex\n",
    "d, p = csgraph.dijkstra(w, indices=[startid], return_predecessors=True, directed=False)\n",
    "print(d.shape)\n",
    "print(p.shape)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 20,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "110"
      ]
     },
     "execution_count": 20,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "np.isinf(d[0]).sum()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 21,
   "metadata": {},
   "outputs": [],
   "source": [
    "path = []\n",
    "node_id = endid\n",
    "while True:\n",
    "    path.append(node_id)\n",
    "    if node_id == startid or node_id < 0:\n",
    "        break\n",
    "    node_id = p[0, node_id]\n",
    "path = np.array(path)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 22,
   "metadata": {},
   "outputs": [
    {
     "ename": "ValueError",
     "evalue": "assignment destination is read-only",
     "output_type": "error",
     "traceback": [
      "\u001b[1;31m---------------------------------------------------------------------------\u001b[0m",
      "\u001b[1;31mValueError\u001b[0m                                Traceback (most recent call last)",
      "\u001b[1;32m<ipython-input-22-8e148d29628c>\u001b[0m in \u001b[0;36m<module>\u001b[1;34m()\u001b[0m\n\u001b[0;32m      1\u001b[0m \u001b[1;31m#%figonly=用dijkstra计算最短路径\u001b[0m\u001b[1;33m\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0;32m      2\u001b[0m \u001b[0mx\u001b[0m\u001b[1;33m,\u001b[0m \u001b[0my\u001b[0m \u001b[1;33m=\u001b[0m \u001b[0mpath\u001b[0m \u001b[1;33m%\u001b[0m \u001b[0mW\u001b[0m\u001b[1;33m,\u001b[0m \u001b[0mpath\u001b[0m \u001b[1;33m//\u001b[0m \u001b[0mW\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[1;32m----> 3\u001b[1;33m \u001b[0mimg\u001b[0m\u001b[1;33m[\u001b[0m\u001b[0my\u001b[0m\u001b[1;33m,\u001b[0m \u001b[0mx\u001b[0m\u001b[1;33m,\u001b[0m \u001b[1;33m:\u001b[0m\u001b[1;33m]\u001b[0m \u001b[1;33m=\u001b[0m \u001b[1;36m0\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0m\u001b[0;32m      4\u001b[0m \u001b[0mfig\u001b[0m\u001b[1;33m,\u001b[0m \u001b[0maxes\u001b[0m \u001b[1;33m=\u001b[0m \u001b[0mpl\u001b[0m\u001b[1;33m.\u001b[0m\u001b[0msubplots\u001b[0m\u001b[1;33m(\u001b[0m\u001b[1;36m1\u001b[0m\u001b[1;33m,\u001b[0m \u001b[1;36m2\u001b[0m\u001b[1;33m,\u001b[0m \u001b[0mfigsize\u001b[0m\u001b[1;33m=\u001b[0m\u001b[1;33m(\u001b[0m\u001b[1;36m16\u001b[0m\u001b[1;33m,\u001b[0m \u001b[1;36m12\u001b[0m\u001b[1;33m)\u001b[0m\u001b[1;33m)\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0;32m      5\u001b[0m \u001b[0maxes\u001b[0m\u001b[1;33m[\u001b[0m\u001b[1;36m0\u001b[0m\u001b[1;33m]\u001b[0m\u001b[1;33m.\u001b[0m\u001b[0mimshow\u001b[0m\u001b[1;33m(\u001b[0m\u001b[0mimg\u001b[0m\u001b[1;33m)\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n",
      "\u001b[1;31mValueError\u001b[0m: assignment destination is read-only"
     ]
    }
   ],
   "source": [
    "#%figonly=用dijkstra计算最短路径\n",
    "x, y = path % W, path // W\n",
    "img[y, x, :] = 0\n",
    "fig, axes = pl.subplots(1, 2, figsize=(16, 12))\n",
    "axes[0].imshow(img)\n",
    "axes[1].imshow(bimg, cmap=\"gray\")\n",
    "for ax in axes:\n",
    "    ax.axis(\"off\")\n",
    "fig.subplots_adjust(0, 0, 1, 1, 0, 0)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "> **SOURCE**\n",
    "\n",
    "> `scpy2.scipy.hrd_solver`使用`csgraph`计算华容道游戏【横刀立马】布局步数最少的解法。"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "#%hide\n",
    "%exec_python -m scpy2.scipy.hrd_solver"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.6.4"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 1
}
